//Construct Intensity Measures
use "$data_in/ED_CPT", clear
replace encdate = dofc(encdate)

forval i = 1/13 {
use patient provnpi encdate cpt_`i' using "$data_in/ED_CPT", clear
qui replace encdate = dofc(encdate)
rename cpt_`i' cpt
qui drop if mi(cpt)
qui save "$data_out/temp_cpt`i'", replace
}

clear
forval i = 1/13 {
append using "$data_out/temp_cpt`i'"
}
drop if cpt == "XXXXX"

format encdate %d
sort patient encdate

save "$data_out/caper_cpt", replace

use "$data_out/caper_cpt", clear
gen calendar_year = year(encdate)
gen temp_month = month(encdate)
gen temp_fy = inrange(temp_month, 10,12)
gen fy = calendar_year + temp_fy
drop temp* 
gen mod = ""

rename (cpt calendar_year) (hcpcs year)
sort hcpcs mod year

merge m:1 hcpcs mod year using "$lookup/HCPCS_RVU_crosswalk.dta", keep(1 3) keepus(total_fac total_fac_trans endo_base mult_proc glob_days surg_bilat surg_asst)

replace total_fac = total_fac_trans if mi(total_fac)
drop total_fac_trans
save "$data_out/caper_sample", replace

tostring endo_base, replace
replace endo_base = hcpcs if inlist(hcpcs, "29805", "29830", "29840", "29860", "29870", "31505", "31525", "31526")
replace endo_base = hcpcs if inlist(hcpcs, "31575", "31622", "43200", "43235", "43260", "44360", "44376", "45300", "45330")
replace endo_base = hcpcs if inlist(hcpcs, "45378", "46600", "49320", "50951", "52000", "52351", "57452", "58555")








* a. Multiple procedure adjustment-

*aii. Applying multiple procedure adjustment - sorting by person, date of service, provider group, mpr, and allowed amount with largest first
	* Endoscopy, no change. All subsequent, subtract endo base payment unless base HCPC billed, then $0 for all subsequent

gen mpr = inlist(mult_proc,1,2,3)
gen mpr_image = inlist(mult_proc,4)
gen endo_base_ind = (endo_base == hcpcs)
gen mpr_end = mult_proc == 3



gsort patient encdate mpr endo_base -total_fac

bys patient encdate mpr endo_base: gen endo_count = _n
bys patient encdate mpr: gen mpr_count = _n
bys patient encdate mpr_image: gen image_count = _n

gen mpr_adjust = 1
replace mpr_adjust = .5 if mult_proc == 1 & mpr_count == 2  
replace mpr_adjust = .25 if mult_proc == 1 & mpr_count > 2  
replace mpr_adjust = .5 if mult_proc == 2 & mpr_count > 1  

gen endo_base_rvu = total_fac if endo_base_ind == 1
bys patient encdate mpr_end endo_base: egen endo_base_ind_fam = max(endo_base_ind)
bys patient encdate mpr_end endo_base: egen endo_base_fam_rvu = max(endo_base_rvu)

gen endo_adjust_rvu = 0
replace endo_adjust_rvu = endo_base_rvu if mpr_end == 1 & endo_count > 1 
replace total_fac = 0 if endo_base_ind_fam == 1 & mpr_end == 1 & endo_count > 1 




	* Identify same procedure billed more than once 
bys patient encdate hcpcs  : gen temp_prov = _n
	
bys patient encdate hcpcs  : gen temp_rvu = 0 if _n > 1 & temp_surg == 1


replace total_fac = temp_rvu if !mi(temp_rvu)


collapse (sum) total_fac, by(patient encdate provnpi) fast

save "$data_out/rvu", replace

//only consider aggregate work RVU's



use "$data_out/caper_cpt", clear
gen calendar_year = year(encdate)
gen temp_month = month(encdate)
gen temp_fy = inrange(temp_month, 10,12)
gen fy = calendar_year + temp_fy
drop temp* 
gen mod = ""

rename (cpt calendar_year) (hcpcs year)
sort hcpcs mod year

merge m:1 hcpcs mod year using "$lookup/HCPCS_RVU_crosswalk.dta", keep(1 3) keepus(wrvu)
collapse (sum) wrvu, by(patient encdate provnpi) fast
drop if mi(provnpi)
save "$data_out/wrvu", replace

//BETOS


use "$data_in/ED_CPT", clear
replace encdate = dofc(encdate)
egen num_hcpcs = rownonmiss(cpt_*), strok

//hcpcs
gen hcpcs = 0
forval i = 4/13 {
replace hcpcs = hcpcs + 1 if inrange(substr(cpt_`i',1,1),"A","Z")
}
gen num_cpt = num_hcpcs - hcpcs

 forval i = 4/13 {

rename cpt_`i' cpt
sort cpt
merge m:1 cpt using "$lookup/BETOS_13", keepus(proc) keep(1 3) nogen
rename cpt cpt_`i'
rename proc B`i'
}

local betos = "B4 B5 B6 B7 B8 B9 B10 B11 B12 B13"
egen B_Tests = anycount(`betos'), v(7)
egen B_Procedures = anycount(`betos'), v(6)
egen B_other = anycount(`betos'), v(5)
egen B_Imaging2 = anycount(`betos'), v(4)
egen B_Unclassified = anycount(`betos'), v(3)
egen B_DME = anycount(`betos'), v(1)



keep patient encdate provnpi num_cpt num_hcpcs B_*
bys patient encdate provnpi (num_hcpcs): keep if _n ==_N

drop if mi(provnpi)

merge 1:1 patient encdate provnpi using "$data_out/wrvu.dta", nogen keep(3)
rename patient pid_pde
merge 1:1 pid_pde encdate provnpi using "$data_out/provider_propensity", keep(2 3)  nogen


cap drop residual 
cap drop mean_residual
gen diagnostics = B_Test +B_Imaging
gen total_cpt = diagnostics + B_Proc
gen log_rvu = log(wrvu)


reghdfe log_rvu i.age_bin race_white female junior_enlisted college married longevity afqt_p, absorb(hymd diag mos) res(rvu_residual)

reghdfe total_cpt i.age_bin race_white female junior_enlisted college married longevity afqt_p, absorb(hymd diag mos) res(total_residual)
reghdfe diagnostics i.age_bin race_white female junior_enlisted college married longevity afqt_p, absorb(hymd diag mos) res(diag_residual)

egen mean_residual_rvu = mean(rvu_residual), by(provID year)

egen mean_residual_diag = mean(diag_residual), by(provID year)
egen mean_residual_total = mean(total_residual), by(provID year)


bys provID year: gen N = _N

gen rvu_propensity = (mean_residual_rvu * N - rvu_residual)/(N-1)

gen diag_propensity = (mean_residual_diag * N - diag_residual)/(N-1)
gen total_propensity = (mean_residual_total * N - total_residual)/(N-1)


save "$data_out/Intensity", replace









use "$data_out/wrvu.dta", clear

drop if mi(provnpi)
rename patient pid_pde
merge m:1 pid_pde encdate provnpi using "$data_out/provider_propensity", keep(2 3)  nogen

recode wrvu .=0

gen log_rvu = log(wrvu)


//BETOS


use "$data_in/ED_CPT", clear
replace encdate = dofc(encdate)
egen num_hcpcs = rownonmiss(cpt_*), strok

//hcpcs
gen hcpcs = 0
forval i = 4/13 {
replace hcpcs = hcpcs + 1 if inrange(substr(cpt_`i',1,1),"A","Z")
}
gen num_cpt = num_hcpcs - hcpcs

 forval i = 4/13 {

rename cpt_`i' cpt
sort cpt
merge m:1 cpt using "$lookup/BETOS_13", keepus(proc) keep(1 3) nogen
rename cpt cpt_`i'
rename proc B`i'
}

local betos = "B4 B5 B6 B7 B8 B9 B10 B11 B12 B13"
egen B_Tests = anycount(`betos'), v(7)
egen B_Procedures = anycount(`betos'), v(6)
egen B_other = anycount(`betos'), v(5)
egen B_Imaging2 = anycount(`betos'), v(4)
egen B_Unclassified = anycount(`betos'), v(3)
egen B_DME = anycount(`betos'), v(1)



keep patient encdate provnpi num_cpt num_hcpcs B_*
bys patient encdate provnpi (num_hcpcs): keep if _n ==_N

drop if mi(provnpi)

merge 1:1 patient encdate provnpi using "$data_out/wrvu.dta", nogen keep(3)
rename patient pid_pde
merge 1:1 pid_pde encdate provnpi using "$data_out/provider_propensity", keep(2 3)  nogen


cap drop residual 
cap drop mean_residual
gen diagnostics = B_Test +B_Imaging
gen total_cpt = diagnostics + B_Proc
gen log_rvu = log(wrvu)


reghdfe log_rvu i.age_bin race_white female junior_enlisted college married longevity afqt_p, absorb(hym hd diag mos) res(rvu_residual)

reghdfe total_cpt i.age_bin race_white female junior_enlisted college married longevity afqt_p, absorb(hym hd diag mos) res(total_residual)
reghdfe diagnostics i.age_bin race_white female junior_enlisted college married longevity afqt_p, absorb(hym hd diag mos) res(diag_residual)

egen mean_residual_rvu = mean(rvu_residual), by(provID year)

egen mean_residual_diag = mean(diag_residual), by(provID year)
egen mean_residual_total = mean(total_residual), by(provID year)


bys provID year: gen N = _N

gen rvu_propensity = (mean_residual_rvu * N - rvu_residual)/(N-1)

gen diag_propensity = (mean_residual_diag * N - diag_residual)/(N-1)
gen total_propensity = (mean_residual_total * N - total_residual)/(N-1)


save "$data_out/Intensity", replace



	foreach outcome in scripts_365 prov7  day_180  {
	eststo `outcome'_C:	ivreghdfe `outcome' $spec1
	eststo `outcome'_diag:	ivreghdfe `outcome' i.age_bin race_white female junior_enlisted college married longevity afqt_p diag_propensity (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) 
	eststo `outcome'_BETOS:	ivreghdfe `outcome' i.age_bin race_white female junior_enlisted college married longevity afqt_p total_propensity (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) 
	eststo `outcome'_RVU:	ivreghdfe `outcome' i.age_bin race_white female junior_enlisted college married longevity afqt_p rvu_propensity (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) 



}

merge 1:1 pid_pde encdate using "$data_out/Promotions.dta", keep(3) nogen
eststo promotion1_C:	ivreghdfe promotion1 $spec_promotions1
eststo Promotion_diag:	ivreghdfe promotion1 diag_propensity $spec_promotions1
eststo Promotion_BETOS:	ivreghdfe promotion1 total_propensity $spec_promotions1
eststo Promotion_RVU:	ivreghdfe promotion1 rvu_propensity $spec_promotions1



 merge 1:1 pid_pde encdate using "$data_out/Flags.dta", keep(3) nogen
 
eststo Flag_C:	ivreghdfe adverse_flag $spec1
eststo Flag_diag:	ivreghdfe adverse_flag diag_propensity $spec1
eststo Flag_BETOS:	ivreghdfe adverse_flag  total_propensity $spec1
eststo Flag_RVU:	ivreghdfe adverse_flag  rvu_propensity $spec1


//Tables
label var ed_prescription "Prescriptions in 1 Year"
	
	#delimit ;

esttab scripts_365_C scripts_365_diag scripts_365_BETOS scripts_365_RVU using "Intensity.rtf" , 
	replace label se star(* 0.10 ** 0.05 *** 0.01)
	b(4) se(4)
	stats( ,
	fmt(a2 %-9.0fc))
	nonotes
		keep(ed_prescription)
		mtitles("Main Results" "Controlling for BETOS Number of Tests and Images Propensity" "Controlling for BETOS Number of Tests, Images, and Procedures Propensity" "Controlling for Log of Work RVU Propensity")
		title(Controlling for Intensity of Medical Treatment During ED Visit);
#delimit cr

label var ed_prescription "7 Providers in 1 Year"
	
	#delimit ;

esttab prov7_C prov7_diag prov7_BETOS prov7_RVU using "Intensity.rtf" , 
	append label se star(* 0.10 ** 0.05 *** 0.01)
	b(4) se(4)
	stats( ,
	fmt(a2 %-9.0fc))
	nonotes
	nomtitles
	nonum
		keep(ed_prescription);
#delimit cr


label var ed_prescription "180 Days of Supply in 1 Year"
	
	#delimit ;

esttab day_180_C day_180_diag day_180_BETOS day_180_RVU using "Intensity.rtf" , 
	append label se star(* 0.10 ** 0.05 *** 0.01)
	b(4) se(4)
	stats( ,
	fmt(a2 %-9.0fc))
	nonotes
	nomtitles
	nonum
		keep(ed_prescription);
#delimit cr

label var ed_prescription "Promotion in 1 Year"
	
	#delimit ;

esttab promotion1_C Promotion_diag Promotion_BETOS Promotion_RVU using "Intensity.rtf" , 
	append label se star(* 0.10 ** 0.05 *** 0.01)
	b(4) se(4)
	stats( ,
	fmt(a2 %-9.0fc))
	nonotes
	nomtitles
	nonum
		keep(ed_prescription);
#delimit cr

label var ed_prescription "Adverse Flag in 1 Year (Army Enlisted Only)"
	
	#delimit ;

esttab Flag_C Flag_diag Flag_BETOS Flag_RVU using "Intensity.rtf" , 
	append label se star(* 0.10 ** 0.05 *** 0.01)
	b(4) se(4)
	stats( ,
	fmt(a2 %-9.0fc))
	nonotes
	nomtitles
	nonum
		keep(ed_prescription);
#delimit cr


